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ABSTRACT 

Assuming the dark matter halo of the Milky Way as a non-spherical potential 
i-C ■ (i.e. triaxial, prolate, oblate), we show how the assembling process of the Milky 

Way halo, may have left long lasting stellar halo kinematic fossils only due to the 
shape of the dark matter halo. In contrast with tidal streams, associated with 
recent satellite accretion events, these stellar kinematic groups will typically show 
inhomogeneous chemical and stellar population properties. However, they may 
be dominated by a single accretion event for certain mass assembling histories. If 
the detection of these peculiar kinematic stellar groups is confirmed, they would 
be the smoking gun for the predicted triaxiality of dark halos in cosmological 
galaxy formation scenarios. 

Subject headings: Galaxy: halo — Galaxy: kinematics and dynamics — Galaxy: 
structure — 



1. Introduction 

The standard cosmological model (ACDM) has reached a development stage in which 
cosmological tests at the scale of galaxies become possible (e.g. Van den Bosch 1998; Courteau 
& Rix 1999; Gnedin et al. 2006; Pizagno et al. 2007). However, a lack of understanding in 
the complicated physics that mediates the evolution of baryonic matter has meant that these 
tests are hard to implement, as their application depends on this physics. In some cases, 
controversial results are found out of these studies (e.g. Klypin et al. 1999; Moore 1994), 
stimulating suggestions of modifications for the model, but probably also the need of more 
accurate comparisons between theory and observations (e.g. Valenzuela et al. 2007; Simon 
& Geha 2007). 
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Under the ACDM picture of structure formation, galactic sized halos are assembled by 
the accretion of smaller structures or subhalos (Kauffmann et al. 1993; Ghigna et al. 1998; 
Klypin et al. 1999), or by major mergers, leaving imprints from these processes in the dark 
matter kinematics and leading to a triaxial halo shape. Besides the assembling process, 
dark matter halos triaxiality may be produced by dynamical secular evolution like the radial 
orbit instability (Aguilar & Merritt 1990; Barnes et al. 2005). All together triaxiality is an 
ubiquitous property of dark matter structures in a hierarchical universe (Allgood et al. 2006; 
Vera-Ciro et al. 2011). Studies aiming to detect triaxiality, have to face the complication 
of halo shape evolution or even that triaxialty may be erased or turned into oblateness by 
the baryonic galaxy formation, as a consequence, a robust detection or rejection of a non- 
spherical halo is still a valuable result constraining the efficiency of dynamical mechanisms 
of halo triaxiality evolution. Different observational strategies have been proposed in order 
to quantify triaxiality in galactic halos. One of them is based on the Sagittarius stellar tidal 
stream morphology and kinematics (Law et al. 2009), other on hypervelocity stars (Gnedin 
et al. 2005), one more on galaxy satellite systems distribution (Zentner et al. 2005). However, 
a clear observational detection of dark matter halo non-sphericity is still a challenge, mostly 
because of the need of higher precision in observations and because of degeneracies triggered 
by the baryonic galaxy structure. Based on the previous discussion our paper is focused on 
detecting a non-spherical halo, therefore we will use triaxiality in a very broad sense without 
distinguishing between true triaxiality, prolateness or oblateness. However at the end of the 
paper we will briefly discuss possible differences. 

The aim of this paper is to introduce a possible new strategy in order to detect triaxiality 
in the Milky Way halo. Our strategy is a generalization of the discussion presented by 
Penarrubia et al. (2009), for dwarf spheroidal satellite galaxies. Peharrubia et al. (2009) 
shows that globular clusters stellar remnants are long lived in the potential of triaxial dwarf 
galaxies, explaining the origin of cold stellar structures observed in UMinor and other Milky 
Way satellites. The kinematic groups correspond to halo version of galactic disks stellar 
moving groups (Antoja et al. 2009), however the mechanism is not identical. In this work we 
propose that if the Milky Way sits in a triaxial dark matter halo, the stellar halo assembling 
history is able to populate the quasi-resonant orbits triggering kinematic stellar groups. In 
order to prove our thesis, we perform numerical simulations of test particles moving inside 
the potential generated by dark matter halos with different shapes, the population of the 
orbital structure is randomly generated. We analyze the resulting kinematic and orbital 
structure in order to decide if the stellar kinematic structure may be evidence for the halo 
triaxiality. 

This letter is organized as follows. In Section 2 the 3-D galactic halo potential used to 
compute orbits is briefly described. In Section 3 we introduce a strategy aimed to efficiently 
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explore the stellar phase space accessible to an hypothetical observer, we also present the 
results of our numerical simulations. Finally, in Section 4 we present a discussion of our 
results and our conclusions. 



2. The Model for the Galactic Halo 



For simplicity our model assumes that the Galaxy is surrounded by a static dark matter 
halo with the density profile proposed by Navarro et al. (1996). It is important to mention 
that in our study we do not pretend to use a Milky Way state of the art model, instead we 
present a proof of concept that a triaxial halo develops and preserves abundance structure 
in the stellar phase space because of the resonant orbital structure. For example, we do 
not include yet any disk effect even that Solar Neighborhood stellar orbits will be clearly 
influenced. Even more, the Galactic dark matter halo shape has been evolving during time, 
however we still use a static model. Our study does not loose generality because, as we show, 
the dominant halo resonant orbital structure is present only inside a non-spherical halo. In 
any case, our results provide a lower limit to the possible galactic triaxiality, which may be 
comparable to the kinematics sampled by an observer further away from the disk or after 
discarding if possible, any disk effect. We will return to this point in the Discussion section. 

Under our assumptions, the gravitational potential generated by the halo takes the form 



$(x,y,z) = 2nGabcp r 2 s I - - = (1) 

' r s + s(r) v /( a 2 + r) ( 6 2 + r )( c 2 + r) ' 



from Peharrubia et al. (2009). Here, the dimensionless quantities a, b and c, are the three 
main axes, p , is the characteristic density of the halo and r s , is the radial scale. We use 
elliptical coordinates, where, 

x 2 \P" z 2 

or + r ¥ + r c z + r 
Halo triaxiality is given then, by a, b and c. 

In the next section, we explore the effect of halo triaxiality in the resulting orbital 
structure. The model density structure is determined by five free parameters, rs, Po and the 
axis ratios a, b and c. Ideally, all can be adjusted using Milky Way observations, however 
this is a not a finished task in Galactic astronomy and we only seek to distinguish between 
orbits in a spherical and triaxial halo. Instead we adopt rs = 8.5 kpc, in order to guarantee 
a flat rotation curve at the solar radius in the absence of a disk, and p — 2.46 M© pc~ 3 , 
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obtained assuming a maximum rotation velocity of 220 km/s. The possible axes ratios are 
even harder to constrain, thus, we decided to explore four cases: One (a=1.47, b=1.22, 
c=0.98), motivated by the halo model that best fits Sagittarius stream (Law et al. 2009), a 
prolate model, an oblate one and a spherical halo model. 

3. Numerical Test Particle Simulations 

We tracked the orbital motion of a number of stars under the influence of the gravita- 
tional potential generated by a steady triaxial dark matter halo as described in the previous 
section. The equations of motion are numerically solved with the Bulirsh-Stoer method 
(Press et al. 1992). We monitored the individual particle energy constant of motion and 
obtain a maximum relative error in our computations of about 10 -9 . We will analyze two 
different set ups that will enable us to demonstrate the power of orbital structure in order 
to identify a non-spherical halo. After 12 Gigayears , which is comparable to the stellar 
halo formation timescale (Kalirai 2012), we analyzed the kinematic projections vx-vz and 
we use the spectral orbital method developed by Carpintero & Aguilar (1998) in order to 
study the orbital type distribution and investigate if the structure is composed of resonant 
orbits. This method is based on the concept of spectral dynamics introduced by Binney & 
Spergel (1984) that uses the Fourier transform of the time series of each coordinate. The 
method distinguishes correctly between regular and irregular orbits, and identifies various 
families of regular orbits (boxes, loops, tubes, boxlets, etc.), it also recognizes the second-rank 
resonances that bifurcate from them. 



3.1. Satellite Infall Initial Conditions 

Our first set of simulations study the infall of a satellite or globular cluster into the 
dark matter halo potential, near to a resonant orbit (2:3). An important characteristic of 
a triaxial potential is that resonant orbits are dominant compared to the situation in a 
spherical potential where there is no angular preference for particles to settle down, and 
the abundance of all type of orbits is not dominated by resonances. This circumstance is 
an important ingredient to determine the halo triaxiality. The aim of this set of controlled 
experiments is to show that in the vicinity of resonances, stellar orbits librate and stay for 
practically infinite time. In our simulation, the satellite is made of 10 4 test particles with 
an isotropic velocity dispersion of 10 km/s, the number of particles has no effect on the 
orbital calculation accuracy, however it may impact the significance of kinematic structures, 
we will comeback to this issue later. Satellite self-gravity is neglected, however including this 
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effect will produce only a delay in the satellite disruption time and little or no effect in the 
abundance of kinematic stellar structure. We let the system evolve for a period comparable 
to recent constraints on the stellar halo formation time (12 Gigayears), under the potential 
generated by the dark halo for the two cases mentioned in Section 2. Figure 1 shows the 
resultant position and velocities for the satellite stars. 

In Figure 1 (upper panels), the particle positions indicate that disrupted stars with the 
right energy and angular momentum librate around a resonant orbit during all the simulation 
time for the triaxial halo. In contrast, for the spherical case it is clear that particles are 
considerably more mixed by the end of the simulation. Figure 1 (lower panels) shows how 
the stars in the vicinity of resonances are prone to form kinematic groups for the triaxial 
case in contrast to what it happens in the spherical halo. This is similar to what has been 
recently claimed by Lux et al. (2012) for a galactic globular cluster. These groups would 
look different from different positions with respect to this fish-like orbit although we have 
focused on the difference between the spherical and the triaxial case. Figure 2 (upper panels) 
shows how velocity space would look for example at the solar neighborhood at the end of 
the simulation, if the Sun were close to the "fish head" orbit (left panels for a spherical 
halo and right panels for a triaxial halo). The observer will find two symmetric arm-like 
kinematic structures close to 100 km/s. We assumed the solar position at (x=-8.5 kpc, y=0, 
z=0), and take as stars in the "solar neighborhood", those found inside a sphere of radius 1 
kpc, centered at the solar position, however we could have placed the observer at a different 
position. In order to illustrate that, Fig. 2 (bottom panels) shows also the velocity space, 
but for stars that are inside a sphere centered at (x=-8.5 kpc, y=0, z=20 kpc). 

For both cases, although at different positions, the velocity space shows two arms that 
correspond to stars that are moving along the resonant orbit in opposite directions. Kine- 
matic groups formed in this manner would provide an evidence of a non-spherical halo. 

3.2. Random Velocity Initial Conditions 

Once we provided a good evidence of kinematic stellar structure triggered by reso- 
nances in a triaxial dark matter halo in a single satellite accretion event, we have to asses 
the probability that a long stellar halo assembly history may populate orbits in the vicin- 
ity of resonances triggering long lasting stellar moving groups. We face three difficulties, 
the number of particles is inversely proportional to noise which could artificially mimic the 
kinematic groups, therefore if we want to have a reliable moving group detection we need to 
have a large number of particles. Secondly, an observer will have good kinematic accuracy 
only inside a limited neighborhood, finally during the galaxy assembly history, satellites 
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and globular clusters have been accreted from many different directions. For that reason 
instead of populating the whole halo with particles, we use the following initial conditions 
set-up. We define a sphere with radius R, around a putative observer, we randomly populate 
the sphere with particles, also with randomly selected velocities between zero and the local 
escape velocity, mimicking the contribution of many accretion events with different orienta- 
tions and energies, etc, but assuring that the corresponding orbit is bounded to the galaxy. 
This is just an strategy to explore the available phase space, is not pretending to be a full 
self-consistent generation of position and velocities. Stars librating around resonances will 
comeback eventually to the observer position therefore we populate orbits that the artificial 
observer will detect, increasing our particle statistics at a low computing price. The process 
produce persisting kinematic groups at the velocities around resonances if the process is 
efficient to explore the available phase space, because irregular and open orbits will almost 
never comeback. In contrast a non-resonant region will present nearly evenly distributed 
types of orbits as a result of phase mixing, and also non-kinematic groups after some mixing 
time scales. 

We set 2 x 10 6 stars randomly distributed inside the "observer neighborhood" , defined 
as mentioned above, we assigned random velocities between zero and the escape velocity. We 
removed from the initial conditions the stars that were moving close to the X-Y plane. These 
stars remain in this plane all along the simulation and will always have vz close to zero. If we 
have populated the whole halo this orbits concentration would be compensated with neighbor 
orbits thrown out from different galaxy positions. Although this orbits would appear in the 
velocity plane as an horizontal band, they do not form a kinematic group, we verified this 
using an spectral orbit clasification (Carpintero & Aguilar 1998), as discussed below. We let 
the system evolve during 12 Gigayears and at the end of this integration time, we studied 
the stellar kinematic distribution. This integration time is long enough to be comparable to 
the Milky Way halo evolution (Kalirai 2012). For our purpose, we consider the stars that 
are found after the simulation, in the "artificial observer neighborhood". Figure 3 shows 
the resulting kinematic distribution of stars after performing the simulation for two different 
cases. Upper panel of Fig. 3 corresponds to a spheric halo (a=1.0, b=1.0, c=1.0) and bottom 
panel corresponds to a triaxial halo (a=1.47, b=1.22, c=0.98). 

From Figure 3, the velocity space for the triaxial halos shows two arms in the same 
position as in the case of the falling satellite. The fact that the kinematic structures do not 
appear in the spheric case strongly suggests that these kinematic arms correspond to stars 
librating in the vicinity of a resonance. In order to proof our interpretation, we use the 
spectral orbital method mentioned in Section 3 and show that the more outstanding kine- 
matic structures have an important contribution of resonant orbits, in contrast a featureless 
region has an almost even distribution of orbital types, and the spherical case shows mostly 
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open orbits. Although we do not pretend to assess yet the detectability of the kinematic 
structure triggered by a non-spherical halo we would like to show that the importance of 
resonant orbits in our velocity structure is not an artifact of small number statistics or poor 
sampling. We would like to stress that the orbital type histogram presented in figure 3 shows 
that particles in the 2:3 resonance are 38.5% of the total, significantly larger than the typ- 
ical fluctuation amplitude. Similar calculation were performed using an oblate and prolate 
models, results although non identical they are qualitatively the same as the ones presented 
in this paper. The corresponding details will be presented in a future study. 

Fig. 4 shows also the velocity space for the triaxial halos, but for stars that are inside a 
sphere centered at a large distance from the galactic disk. Notice that the arm-like features 
are still detected, just as in the satellites case. Therefore, details of the observers position 
and galactic model, through important are not key in order to detect the kinematic groups. 

4. Discussion and Conclusions 

Using test particle simulations, we studied the motion of stars under the influence of 
the gravitational potential generated by a steady triaxial dark matter halo and a spherical 
one. Our simulations show that a non-spherical shape of the Milky Way dark matter (tri- 
axial, prolate, oblate) halo determines an important amount of the stellar halo kinematic 
structure, and that this is a long lasting feature. There is some dependence on the kinematic 
structure distinguishing the triaxial case from the prolate and oblate situation, however the 
most outstanding difference is the comparison with the spherical case, where no kinematic 
structure is triggered. 

In a satellite accretion event into the triaxial halo potential, close to a resonant orbit, star 
particles librate around resonances for times longer than the Universe age (16 Gyrs). Even if 
the stars have initial random velocities, mimicking random satellite accretion events during 
the Galactic halo assembling, many of the stars still stay around resonant orbits and are prone 
to form kinematic groups if the dark matter halo is triaxial. The reason is that an important 
fraction of stars will be librating around resonances and they will be prone to form kinematic 
groups (See Figs. 3 and 4). These groups might be observationally detectable, however this 
endeavor will face important challenges. First, the disk contribution to the gravitational 
potential may trigger orbital structure as it has been discussed by Lux et al. (2012), however 
either the disk effect may be characterized and avoided, or the search may focuse on stars 
further away from the disk effect in order to minimize the disk graviational effect. The era 
of GAIA will provide us with strong constraints to the Galaxy stellar structure, as well as 
6D spatial information and chemical composition of a huge number of stars and at distances 
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further away from the solar neighborhood. Careful analysis of these data compared with 
the results of our simulations may shed some light on the shape of the Milky Way halo. A 
single event such as the tidal disruption of a globular cluster or a satellite, may produce long 
lasting kinematic groups even in configuration space (Lux et al. 2012; Penarrubia et al. 2009) 
with some degeneracy in their interpretation due to the dark matter halo and the dynamical 
age of the accretion event. However, we will be able to distinguish these groups from those 
associated with the dark matter halo shape with the help of the stellar population data, 
which must be very homogeneous. The detection of these type of kinematic groups, would 
be a clear indication of the triaxiality (non-sphericity) of the Milky Way dark matter halo, 
and will be of great importance for the galaxy formation theories and the ACDM scenario. 

We thank Antonio Peimbert for some enlightening discussions. We thank DGAPA- 
PAPIIT through grants IN110711, IN117111-3 and CONACyT grant CB-128556. 
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Fig. 1. — Position and velocities of the satellite stars after 16 Gigayears for the spherical 
(left) and the triaxial halo (right). It is remarkable the small mixing shown after the long 
integration time for particles close to the 2:3 resonant orbit in the triaxial case. 
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Fig. 2. — Velocity space in the solar neighborhood (upper panels) for both cases spherical 
(left panel), triaxial (right panel), setting the Sun close to the fish head orbit. Bottom panels 
show the case for an observer further away by 20 kpc for both cases spherical (left panel), 
triaxial (right panel). 
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Fig. 3. — Kinematic structure in spherical (top panel) and triaxial halo (bottom panel) 
models. We show the vx-vz projection of the velocity space measured by a hypothetical 
observer inside the spherical model, the projection is basically featureless, the inset shows 
a histogram of the orbital types, the contribution of resonant orbits is negligible. Right 
hand panel shows the vx-vz diagram corresponding to an observer inside the triaxial halo 
model. We recover the structure related to the fish-like orbits as in the single satellite 
experiment, and other structures, regardless that in this case we randomly populated the 
orbital distribution. We present three insets in the triaxial halo case: the upper right one 
that shows almost a nearly even distribution of all orbital species, the inset in the middle 
shows a distribution clearly dominated by resonant orbits (2:3), the lower left inset is centered 
in another kinematic structure, again the 5:7 resonant is dominant. We conclude that the 
orbital type distribution supports the resonant origin of the kinematic structure and its 
capability as a triaxiality diagnostic. Both systems were evolved for 2 x 10 9 years. 
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Fig. 4. — Same as Figure 3 c) and d), but for a different position (8.5,0,1.5). Regardless of 
the observer height with respect to the disk plane, it is possible to detect the signature of 
the non spherical halo. 



